#include<iostream>
#include <xmmintrin.h> //SSE
#include <emmintrin.h> //SSE2
#include <pmmintrin.h> //SSE3
#include <tmmintrin.h> //SSSE3
#include <smmintrin.h> //SSE4.1
#include <nmmintrin.h> //SSSE4.2
#include <immintrin.h> //AVX、AVX2
#include<pthread.h>//导入成功
#include<windows.h>
#include<stdlib.h>
#include<cstdlib>
#include<malloc.h>
#include<math.h>
using namespace std;
int N = 2000;
//int num;
float** matrix;
const int thread_num = 8;//常量,开八个线程/还尝试了4个线程

//生成不会出错的测试用例
void m_reset(float** m, int n)
{
    for (int i = 0; i < n; i++)
    {
        // 下三角全部设置为0
        for (int j = 0; j < i; j++)
            m[i][j] = 0;
        //对角线设置为1.0
        m[i][i] = 1.0;
        //上三角设置为随机数
        for (int j = i + 1; j < n; j++)
            m[i][j] = rand() % 100;
    }

    //每一行加上比自己下标小的行
    for (int k = 0; k < n; k++)
        for (int i = k + 1; i < n; i++)
            for (int j = 0; j < n; j++)
                m[i][j] += m[k][j];
}

//定义县城数据结构
typedef struct {
    int k;//消去的轮次
    int t_id;//线程id
    int num;//矩阵的阶数
}threadParam_t;

//线程函数做消去操作
void *threadFunc(void* param)
{
    threadParam_t* p = (threadParam_t*)param;
    int k = p->k;//消去的轮次
    int t_id = p->t_id;//线程编号
    int num = p->num;
    for (int i = k + t_id + 1; i < num; i=i+thread_num)
    {
        for (int j = i + 1; j < num; j++)
            matrix[i][j] = matrix[i][j] - matrix[i][k] * matrix[k][j];
        matrix[i][k] = 0;
    }
    pthread_exit(NULL);
    return NULL;
}

void solve(int n)
{
    for (int k = 0; k < n; ++k)
    {
        //主线程做的事情，处理本行数据
        for (int j = k + 1; j < n; j++)
        {
            matrix[k][j] = matrix[k][j] / matrix[k][k];
        }
        matrix[k][k] = 1.0;
        //int worker_count = n - 1 - k;//工作线程数量

        //其他线程做的事，下面行的消去不能有多少行就开多少个，做不到/代价太高，
        pthread_t handles[thread_num];// 创建对应的 Handle
        threadParam_t param[thread_num];// 创建对应的线程数据结构
        //开八个线程
        for (int t_id = 0; t_id < thread_num; t_id++)
        {
            param[t_id].k = k;
            param[t_id].t_id = t_id;
            param[t_id].num = n;
        }

        //开启线程
        for(int t_id = 0; t_id < thread_num; t_id++)
            pthread_create(&handles[t_id],NULL,threadFunc,(void*)(param+t_id));

        //结束线程
        for(int t_id = 0; t_id < thread_num; t_id++)
            pthread_join(handles[t_id],NULL);

    }


}


int main()
{
    matrix = new float* [N];
    for (int i = 0; i < N; i++)
    {
        matrix[i] = new float[N];
    }

    int size = 32;//数据集增长速度
    long long head1, tail1, freq1;//记录串行时间
    for (int n = size; n < N; n += size)
    {
        m_reset(matrix, n);

        //串行
        QueryPerformanceFrequency((LARGE_INTEGER*)&freq1);
        QueryPerformanceCounter((LARGE_INTEGER*)&head1);
        QueryPerformanceCounter((LARGE_INTEGER*)&tail1);
        int count1 = 0;
        while ((tail1 - head1) < 0.1 * freq1) {
            solve(n);
            QueryPerformanceCounter((LARGE_INTEGER*)&tail1);
            count1++;
        }
        float time1 = (tail1 - head1) * 1000.0 / (freq1 * count1);

        cout << "问题规模 " << n << " 动态平均时间 " << time1 << " ms";
        cout << endl;
        if (n > 320) size = 128;
    }
}